First Order Superconductor to Insulator Transition: Evidence for a Supersolid Phase 
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The superconductor to insulator transition in the presence of long-range Coulomb interactions is 
studied using Monte Carlo techniques. At integer filling the transition is second order. At half 
filling we find evidence of a first order transition to a state with non-zero S',r,ir. No intervening 
supersolid phase is observed. Off half-filling two separate transitions are found. The transition 
from the superconducting to a supersolid phase is second order followed by another second order 
transition between the supersolid and an insulating phase. 
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If a phase possesses both diagonal long range order 
(DLRO), associated with a structural ordering, as well 
as off diagonal long range order (ODLRO), indicating 
phase coherence, it is xeferred to as a supersolid. Initially 
Penrose and OnsagerEl showed that such a phase cannot 
occur if each lattice site is always occupied. However, 
it was later realized that if the occupancy is not an in- 
teger, and if the wave-function is properly symmetrized, 
a supersolid phase can exista. Analogous phases m spin 
models were investigated and shown to be possiblcfl. Ex- 
perimental evidence far a supersolid phase in solid *He 
seems at present scanta although recent experiments have 
seen indications of jiovel behaviorEl possibly relatecLto a 
superhexatic phaseQ. Exploiting a duality argunientl3 the 
analogue of a supersolid phase of vorlex lines in type- 
II superconductors has been proposedQ. Recent numer- 
ical studies have found ^evidence of supersolid phases in 
quantum lattice modelsLTEJ. These studies indicate that 
interactions longer range than on-site are necessary to 
stabilize a supersolid phase. However, true long-range 
interactions were not considered, and the order of the 
transition was not investigated numerically. 

In this communication we study the phase diagram of 
the 2D charged lattice Bose gas at several different fill- 
ing factors. This model should be applicable to short 
coherence length superconductors where the individual 
Cooper PpifS have formed prior to the superconducting 
transitior£3. At integer fillings we find a second order 
transition with SELX K-like exponents in agreement with 
predictions by Mali3. At half-filling we find a single first 
order transition from the superconducting phase into an 
insulating phase with checkerboard charge ordering. No 
supersolid phase is found. Slightly below half-filling wc 
find a second order transition from the superconducting 
phase into a supersolid phase with ODLRO and charge 
ordering followed by another second order transition into 
an insulating phase. A first order supersolid transition 
has recently been observed in a soft sphere modeO. A 
detailed account of our results will be published else- 
wherecJ. ■— . 

Our starting point is the boson Hubbard modetj 
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Here U is the on-site repulsion, fi is the chemical poten- 
tial, z the number of nearest neighbors, and fir — ^l^r 
is the number operator on site r, with <I>J, boson cre- 
ation and annihilation operators. The hopping strength 
is given by t and G is the 1/r Coulomb potential modi- 
fied for the periodic lattice. We shall use either an Ewald 
or lattice Greens function for G with e*^ = l/f^ 
model can be mapped onto a link-current modelT 
defined in terms of integer currents J — [ J^ , , J'^) on 
the links of the (2-|-l)D lattice where r = {x^y) and t is 
the imaginary time which has been discretized. Note that 
J'^ is just the particle density. The integer-link variables 
(J^, jy , J'^) take integer values from —oo to oo, subject 
to the divergence-free constraint given by the continuity 
equation, V • J = 0. If we include long range Coulomb 
interactions we caiuwrite down a purely real action in 
this representation^ 

^ = :^E{^J'('-'^)-^'^'(^'^)} 

r,r 

+ ^Y.^r^v,r)-n,]G{v-v')[r{v',T)^n,]. (2) 

r,r' ,T 

The calculation is performed in the canonical ensemble. 
Hence, the term involving /i is a constant and is only in- 
cluded for computational convenience. At integer filling 
in the absence of long-range interactions this model is in 
the same universality clas&-as_the 3D XY model and has 
an inverted XY transitiontj'Ej. Including long-range in- 
teractions previous studies have fomd that the transition 
can be either first or second ordered at integer filling. 

Having arrived at an equivalent classical model we can 
employ standard Monte Carlo methods to study the in- 
sulating transition by letting K take the place of a "tem- 
perature" . See Ref. O for a detailed description. We 
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first locate the approximate position of the critical point, 
Kc, with low statistiGS-|Calculations. Then the histogram 
Monte Carlo methodcl is used to obtain high precision 
results for a range of K close to Kc. For each system size 
3 - 10 X 10"* Monte Carlo sweeps (MCS) of the entire lat- 
tice are discarded for equilibration followed by 3— 10 x 10"* 
MCS for calculations with a measurement every 10 MCS. 
In order to calculate the error-bars we repeat the entire 
calculation 50 times for each lattice size. Two replicas 
with different random number sequences are run in par- 
allel allowing overlaps between them to be calculated. 
Techniques used in the study of spin glasses then allow 
us to, check that we are calculating equilibrium quanti- 
tiesB. Thus we use in total 6 - 20 x 10^ MCS to build the 
final probability distribution for use in the histograms. 

In order to study the structural ordering we define the 
structure factor, S'k, as follows: 
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Note that S'k is positive definite. At half filling we expect 
the dominant ordering to be associated with the wave- 
vector (tt, tt) corresponding to a checkerboard ordering 
with basis vectors (a, a) and (a, —a). The onset of the 
supercopducting phase will be associated with a non-zero 
stiffnessEj defined as 
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Note tiiat this stiffness corresponds to a non-zero Drude 
weighted. Invoking finite-size scaling we arrive at the 
following finite-size scaling form for the stiffness 
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Here z is the dynamical critical exponent defined through 
the relation ^t- ~ and 5 = K—Kc is the distance to the 
critical point, Kc- Assuming that the transition is second 
order, we see From Eq. (||) that if the aspect ratio, c = 
Lt/L^, is kept constant L'^^^^^p should be independent 
of L at the critical point. In the-jjresence of long-range 
Coulomb interactions Fisher et alnJ has argued that z = 1 
and we shall in the following take that as our starting 
assumption. 

Integer filling. For a filled lattice we do not expect 
any non-trivial charge ordering and there should only 
be a single transition between the superconductor and a 
Mott insulating^hase, governed by a dynamical critical 
exponent z = 1l3. We therefore locate the critical point 
by calculating Lp using cubic samjples. A well defined 
crossing is found at Kc — 0.289(2)113. Using the finite 
size scaling form Lp — p{SL^^'^) a scaling plot of Lp can 
be used to determine v. We find v ~ 0.67113. We also 
calculate the correlation functions both in the time and 



space direction. This yields an independent calculation 
of the dynaroical critical exponent, z, and the anomalous 
dimension if^. By fitting the correlation functions at 
the critical point to power law forms we find z = 1.0 ± 
0.1, rj = —0.15 in reasonable agreement with the 

prediction by Mall3 that the exponents should remain 3D- 
XY like due to the irrelevance of Coulomb interactions 
for this filling factor. However, critical amplitudes could 
possibly be affected by the Coulomb interactions and a 
priori we cannot exoect to find the same value of the 
universal resistivitynJ a± the critical point as was found 
for the 3D-XY modeltZl. If we calculate the frequency 
dependent stiffness: 
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the conductivity can be calculated from the relation 
a{iujn) — 2'KaQp{iuJn)/wn- Here the quantum of con- 
ductance is ctq = {2e)^/h = Rq and LUn = 21m/ Lr 
is the nth Matsubara frequency. With a 1/L finite size 
correction it is possible to collapse all the data for dif- 
ferent lattice sizes onto a single curve at Kc and we find 
R*^/Rq = 1.5 ±0.4. 

Half-integer filling. We now turn to a discussion of 
the half filled case. In this case we presumably also have 
z — lollll, however, none of the other exponents are 
known. In order to locate the structural transition as- 
sociated with the assumed (tt, tt) ordering we calculate 
a cumulant ratio of the structure factor S.„_-„ defined as: 
g{K) = (3- < Sl^^ > / < S^,^ >2)/2. With this defini- 
tion g{K) should approach 1 in the ordered phase and 
in the disordered phase. At the critical point g{K) should 
be independent of L and curves representing g{K) for dif- 
ferent lattice sizes should therefore intersect close to Kc 
independent of the order of the transitioncJ. Our results 
are shown in Fig. |^ for the lattice sizes L = 6, 8, 10, 12, 14. 
The curves for different L intersect at Kc = 0.2559(4). 
However, above Kc the curves for the larger L devel- 
ops a minimum and starts to increase. This behav- 
ior is-usually associated with a first order phase tran- 
sitions. A plot of the probability distribution, P{E) for 
the larger lattices shows a double peak structure-dose 
to Kc = 0.2559(4). Following Lee and KosterlitzH we 
consider A{E,L) = -\nP{E,L). Here P{E,L) is cal- 
culated at Kc{L), the coupling at which the probability 
distribution has two peaks of equal heights at E = Ei, E2 
surrounding a minimum at Em- One can then define a 
bulk free-energy barrier between the ordered and disor- 
dered states as AF{L) = A{Em,L) - A{E„ L). In Fig. | 
we show results for A{E) for L = 8, 10, 12, 14. As a 
function of L, AF(L) increases monotonically. However, 
since AF{L) haSr^ot approached the asymptotic form 
AF{L) ^ i^d+z-v^ expect that we in all cases have 
L < ^. S'tt 7r is significantly different in the two minima. 

In Fig. |l| we show results for the cumulant ratio asso- 
ciated with the winding number and hence the onset of 
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superfluidity. This quantity as well as Lp shows that the 
superfluid order parameter disappears at the same Kc 
where the first order structural transition takes place. 
We therefore conclude that the action Eq. (||) has a sin- 
gle first order phase transition at half-filling. We believe 
that this implies that the underlying quantum Hamilto- 
nian therefore also has a first-order transition, although 
it is possible to argue that neglected irrelevant terms in 
Eq. (|) could cause the transition to be second order 
for the quantum problem. However, we have explicitly 
checked that the transition remains first order for both 
the Ewald and lattice Greens function form of the poten- 
tial. The above scenaro is in good agreement with recent 
work by Batrouni et al., and Scalettar et al. [Tl| . 

Off half-filling. As the range of the interaction is in- 
creased a particulary rich phase diagram emerges for the 
boson Hubbard modeEj. At half-filling the insulating 
phase at small K corresponds to a checkerboard charge 
ordering with unit cell twice as large as the underlying 
lattice. Due to the limited lattice sizes that it is pos- 
sible to study numerically, the most convenient filling 
factor away from half-filling is 7/16. At this filling factor 
we find, in the limit ii' — > 0, an insulating phase with 
charge ordering corresponding to vacancies in a (0,4a), 
(4a, 2a), (4a, —2a) triangular lattice superimposed on the 
half-filled (tt, tt) orderingij. Analogous orderings can be 
found at other filling factors. This charge ordering is 
commensurate with the lattice sizes L = 8, 16. Thus, 
at this filling factor we can use the S^^tt structure factor 
to monitor any charge ordering even though the actual 
ground-state corresponds to xtrdering at more than one 
wave-vector. Scalettar et alllil have shown that z — 1 
at the superfluid-supersolid transition. Due to the long 
range interactions we expect z = 1 also at the supersolid- 
solid transition as opposed to the short-range caseEil. In 
our simulations we therefore use z — 1 throughout. Our 
results are shown in Fig. |[ Two separate transitions 
are visible at = 0.2125(15) and Kc2 = 0.2305(15) 
surrounding a supersolid phase. At the transition be- 
tween the superconducting and supersolid phase we find 
l/jy ~ 1.8±.3. The transition between the supersolid and 
insulating phase is tentatively characterized by the expo- 
nent l/v ^ 1.1 ± 0.3, but a detailed analysis is hindered 
by the fact that the L = 12 system is not commensu- 
rate with the assumed ground-state. We have also calcu- 
lated the structure factors 5'7r,o and Sq^t^ corresponding 
to coUinear ordering. In none of the above cases have we 
observed charge ordering at these wave vectors. 

In conclusion, we have found a /jrst- order superconduc- 
tor to insulator transition at half-filling. Off half filling 
we find evidence for a supersolid phase surrounded by 
two second order phase transitions. 
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FIG. 1. The cumulant ratio of the structure factor, Stt.tt, 
and winding number, W , as a function of K at half-filhng. 
The lines indicate results obtained using the histogram Monte 
Carlo method for different lattice sizes. The circles, squares 
and triangles denote seperate calculations with lesser preci- 
sion for lattices of size L = 6, 8, 10, respectively. 
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FIG. 2. A{E) = - \nP{E) as a function of E/L'-^ for the 
lattice sizes L = 8, 10, 12, 14 at half-filling. P{E) was cal- 
culated at Kc{L) = 0.25345, 0.25465, 0.25504, .25548, respec- 
tively. The relative height of the central peak increases with 
L. 
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FIG. 3. Results for a filling factor of 7/16. Shown are the 
the cumulant ratio for the structure factor, and Lp as a func- 
tion of K. 
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